*-----------------------------------------------*
* Title: Replication Code for Jeannet, A.-M. and Dražanová, L. 2023. Cohort Differences in Attitudes Towards Sexual Orientation: The Formative Political Climate as a Socializing Agent. Frontiers in Political Science 5. doi: 10.3389/fpos.2023.1223274   *
* Author: Lenka Dražanová                            *
* Date: 03/10/2023                                   *
* Description: This .do file provides all the necessary steps to replicate the analysis in the article  *
*-----------------------------------------------*

*use the European Social Survey's integrated file for 10 waves (2002-2020) "ESS10 - integrated file" 10.21338/ess10e03_1 available at https://ess-search.nsd.no/en/study/172ac431-2a06-41df-9dab-c1fd8f3877e7

local datafile "data/ESS1-10.dta"
use `datafile', clear


** Section 1: Data Preparation

**COUNTRY SUBSAMPLE
encode cntry, gen (country)
***keeping only countries needed for the analysis: Austria, Belgium, Denmark, Finland, France, Germany, Ireland, Italy, Netherlands, Norway, United Kingdom, Sweden, Switzerland
egen OK = anymatch(country), values(2 3 5 8 9 12 13 14 18 21 27 28 34)
keep if OK
*separating the two Germanies
recode country 8=40 if intewde==1
label define country 40 "DE-E", modify
drop if country==40

*DEPENDENT VARIABLE (Gays and lesbians free to live life as they wish)
***recoding the dependent variable so that higher numbers mean more positive attitudes
recode freehms 1=5 2=4 3=3 4=2 5=1, gen(freelife)

*INDIVIDUAL LEVEL INDEPENDENT VARIABLES
**AGE AND COHORT VARIABLES
***creating variable allage
gen allage = .m
replace allage=agea
***making sure there are enough observations for all ages
drop if allage<18
drop if allage>85 & allage<.
gen agesq=allage*allage
***Creating cohorts
generate cohort=.m
replace cohort=1 if yrbrn<1926
replace cohort=2 if yrbrn>=1926 & yrbrn<=1930
replace cohort=3 if yrbrn>=1931 & yrbrn<=1935
replace cohort=4 if yrbrn>=1936 & yrbrn<=1940
replace cohort=5 if yrbrn>=1941 & yrbrn<=1945
replace cohort=6 if yrbrn>=1946 & yrbrn<=1950
replace cohort=7 if yrbrn>=1951 & yrbrn<=1955
replace cohort=8 if yrbrn>=1956 & yrbrn<=1960
replace cohort=9 if yrbrn>=1961 & yrbrn<=1965
replace cohort=10 if yrbrn>=1966 & yrbrn<=1970
replace cohort=11 if yrbrn>=1971 & yrbrn<=1975
replace cohort=12 if yrbrn>=1976 & yrbrn<=1980
replace cohort=13 if yrbrn>=1981 & yrbrn<=1985
replace cohort=14 if yrbrn>=1986 & yrbrn<=1990
replace cohort=15 if yrbrn>=1991 & yrbrn<=1995
replace cohort=16 if yrbrn>=1996 & yrbrn<=2000
replace cohort=17 if yrbrn>=2001 & yrbrn<=.

label define cohort 1 "<1926 prewar" 2 " 26-30 post-war" 3 "31-35" 4 "36-40" ///
5 "41-45" 6 "46-50" 7 "51-55" 8 "56-60" 9 "61-65" 10 "66-70" 11 "71-75" 12 "76-80" ///
13 "81-85" 14 "86-90" 15 "91-95" 16 "96-00" 17 ">2001"
label variable cohort cohort

egen miss = rowmiss (cohort)
drop if miss==1
drop miss

***some cohorts have to be excluded from the analysis because they are not enough observations per period
drop if cohort==1
drop if cohort==2
drop if cohort==17

**EDUCATION
generate edu=.m
replace edu=eisced
recode  edu 55=.
replace edu = 1 if eisced== 0 & edulvla == 1 
replace edu = 2 if eisced== 0 & edulvla == 2 
replace edu = 3 if eisced== 0 & edulvla == 3 
replace edu = 5 if eisced== 0 & edulvla == 4 
replace edu = 6 if eisced== 0 & edulvla == 5 
recode edu 0=.
***creating education dummies
recode edu 1=1 .=. *=0, gen(llsec)
recode edu 2=1 .=. *=0, gen(lsec)
recode edu 3=1 4=1 .=. *=0, gen(usec)
recode edu 5=1 .=. *=0, gen(subdeg)
recode edu 6=1 7=1 .=. *=0, gen(uni)


**GENDER
recode gndr 1=0 2=1, gen(female)
**URBAN
recode domicil 1=1 2=1 .=. *=0, gen(urban)
**ECONOMIC SATISFACTION/INCOME DIFFICULTIES
recode hincfel 1=0 2=0 3=1 4=1, gen(incdiff)
**UNEMPLOYED
generate unempl=0
replace unempl= 1 if uempla==1
replace unempl= 1 if uempli==1 
recode uempnyr 1=1 2=1 3=0 4=0 .=., gen(lklunmp)
replace unempl= 1 if lklunmp==1 


*Note: The models will also be using the individual-level independent variables: lrscale and do not need any recoding rlgdgr



*****COHORT LEVEL INDEPENDENT VARIABLES
generate country_num=country
*PRINCIPLES
egen party_id=concat(country_num yrbrn)
merge m:m party_id using "data/cohortprinciples.dta"
drop _merge
**EQUALITY
bysort cohort country: egen meanequalityc = mean(equality)
**TRADITIONALISM
bysort cohort country: egen meantraditionalismc = mean(traditionalism)

*LAWS ABOUT REGISTERED PARTNERSHIPS OR MARRIAGE INDEX
egen gaymlaw_id=concat(country_num yrbrn)

merge m:m gaymlaw_id using "data/gaymlaw.dta"
drop _merge

bysort cohort country: egen meangaymlawc = mean(gaymlaw)

*EDUCATION
bysort cohort country: egen eduunic = mean(100 * (uni==1))

*HISTORICAL UNEMPLOYMENT
egen histunemp_id=concat(country_num yrbrn)

merge m:m histunemp_id using "data/histunempl.dta"
drop _merge

bysort cohort country: egen meanhistunemplc = mean(histunempl)

**NOTE: THE OLDEST COHORT IN EACH COUNTRY ARE ASSIGNED THE VALUE FOR THE YEAR 1956 (THE EARLIEST AVAILABLE YEAR) 
replace meanhistunemplc = 2.893983 if country==2 & cohort==3
replace meanhistunemplc = 2.168675 if country==6 & cohort==3
replace meanhistunemplc = 3.44451 if country==4 & cohort==3
replace meanhistunemplc = 5.343512 if country==10 & cohort==3
replace meanhistunemplc = 0.5896806 if country==13 & cohort==3
replace meanhistunemplc = 0.9915014 if country==14 & cohort==3
replace meanhistunemplc = 1.740644 if country==15 & cohort==3
replace meanhistunemplc = 0.1164325 if country==3 & cohort==3
replace meanhistunemplc = 0.9202191 if country==8 & cohort==3
replace meanhistunemplc = 3.560371 if country==1 & cohort==3
replace meanhistunemplc = 5.563986 if country==5 & cohort==3
replace meanhistunemplc = 1.143678 if country==7 & cohort==3
replace meanhistunemplc = 1.440922 if country==11 & cohort==3
replace meanhistunemplc = 9.90991 if country==12 & cohort==3


*****CREATING PERIOD LEVEL VARIABLES
*PERIOD EDUCATION BY COUNTRY
bysort essround country: egen eduunip = mean(100 * (uni==1))

*PERIOD UNEMPLOYMENT
egen perunempl_id=concat(country_num essround)

merge m:m perunempl_id using "data/perunempl.dta"
drop _merge

*PERIOD GAY MARRIAGE LAWS
egen pergaymlaw_id=concat(country_num essround)

merge m:m pegaymlaw_id using "data/pgaymlaw.dta"
drop _merge


*PRINCIPLES
egen perparty_id=concat(country_num essround)

merge m:m perparty_id using "data/periodprinciples.dta"
drop _merge



* Centering all indivdiual-level continuous independent variables around their grand mean to use in hierarchical models:
foreach x of varlist allage lrscale rlgdgr /// 
{
sum `x', d
gen `x'_cent=`x'-r(mean)
}


*generating missing variable count to use only cases with valid observations across all variables
*capture drop missing 
egen missing=rmiss(allage uni female urban incdiff lrscale rlgdgr ///
meanequalityc meantraditionalismc eduunic meanhomolawc  ///
pequality ptraditionalism eduunip phomolaw ///
)
drop if missing>0

* Section 2: Data Analysis
*we start with a null model

*MODEL 0
mixed freelife || country: || _all:R.cohort || essround: , variance
estimate store M0


outreg2 using cohort_m0.doc, sideway /// 
replace ctitle(Coef) ///
alpha(0.001, 0.01, 0.05) symbol(***, **, *) 


********PREDICTING COHORTS AND PERIOD RANDOM EFFECTS******
*First, we use the egen command with the tag() function to generate three binary 
*indicator variables to tag (set equal to 1) one observation per cohort, one observation
*per period, and one observation per country, respectively
egen pickone_country = tag(country)
egen pickone_cohort = tag(cohort)
egen pickone_p = tag(essround)

	
*COHORTS random effects
egen group = group(country)
gen c02 = .
gen c0se2 = .
gen pr_c_low2 = .
gen pr_c_high2 = .
su group, meanonly

forval g = 1/`r(max)' {
    quietly mixed freelife if group == `g' || _all: R.cohort || essround:, variance
    predict re, reffects relevel(_all)
	predict se, reses relevel(_all)
    replace c02 = re if group == `g' 
	replace c0se2 = se if group == `g' 
    drop re
	drop se
	replace pr_c_low2 = (c02 - c0se2) if group == `g' 
    replace pr_c_high2 = (c02 + c0se2) if group == `g' 
}

********FIGURE 1*********************	
twoway (rarea pr_c_low2 pr_c_high2 cohort, sort fcolor(gs13) lcolor(white)) ///
	(line c02 cohort, sort mcolor(black) msymbol(circle) lcolor(black) ///
	lwidth(medthick)) ///
	, ytitle("") yline(0) ///
	xtitle("") xlabel(3(1)16, labsize(large)) ///
	legend(off) graphregion(fcolor(white) lcolor(white)) ///
	by(country, title("Cohort Random Effects by Country", size(medsmall)))
***************************************	



*PERIOD random effects
gen p02 = .
gen p0se2 = .
gen pr_p_low2 = .
gen pr_p_high2 = .
su group, meanonly

forval g = 1/`r(max)' {
    quietly mixed freelife if group == `g' || _all: R.cohort || essround:, variance
    predict re, reffects relevel(essround)
	predict se, reses relevel(essround)
    replace p02 = re if group == `g' 
	replace p0se2 = se if group == `g' 
    drop re
	drop se
	replace pr_p_low2 = (p02 - p0se2) if group == `g' 
    replace pr_p_high2 = (p02 + p0se2) if group == `g' 
}

********FIGURE 2*********************
twoway (rarea pr_p_low2 pr_p_high2 essround, sort fcolor(gs13) lcolor(white)) ///
	(line p02 essround, sort mcolor(black) msymbol(circle) lcolor(black) ///
	lwidth(medthick)) ///
	, ytitle("") yline(0) ///
	xtitle("") xlabel(1(1)10, labsize(large)) ///
	legend(off) graphregion(fcolor(white) lcolor(white)) ///
	by(country, title("Period Random Effects by Country", size(medsmall)))
***************************************	
	
*MODEL 1 (ONLY INDIVIDUAL LEVEL)
mixed freelife allage_cent i.uni i.female i.urban i.incdiff lrscale_cent rlgdgr_cent ///
|| country: || _all:R.cohort || essround: , variance
estimate store M1

* Save regression results as a table in a .doc file
outreg2 using cohort_m1.doc, sideway /// 
replace ctitle(Coef) ///
alpha(0.001, 0.01, 0.05) symbol(***, **, *) 


*MODEL 2
mixed freelife allage_cent i.uni i.female i.urban i.incdiff lrscale_cent rlgdgr_cent ///
meanequalityc ///
pequality ptraditionalism eduunip perunemp pgaymlaw ///
|| country: || _all:R.cohort || essround: , variance
estimate store M2

* Save regression results as a table in a .doc file
outreg2 using cohort_equalityp.doc, sideway /// 
replace ctitle(Coef) ///
alpha(0.001, 0.01, 0.05, 0.1) symbol(***, **, *, +)


*MODEL 3
mixed freelife allage_cent i.uni i.female i.urban i.incdiff lrscale_cent rlgdgr_cent ///
meantraditionalismc ///
pequality ptraditionalism eduunip perunemp pgaymlaw ///
|| country: || _all:R.cohort || essround: , variance
estimate store M3

* Save regression results as a table in a .doc file
outreg2 using cohort_traditionalismp.doc, sideway /// 
replace ctitle(Coef) ///
alpha(0.001, 0.01, 0.05, 0.1) symbol(***, **, *, +)


*MODEL 4
mixed freelife allage_cent i.uni i.female i.urban i.incdiff lrscale_cent rlgdgr_cent ///
meanequalityc meantraditionalismc ///
eduunic meanhistunemplc meangaymlawc ///
pequality ptraditionalism eduunip perunemp pgaymlaw ///
|| country: || _all:R.cohort || essround: , variance
estimate store M4

* Save regression results as a table in a .doc file
outreg2 using cohort_allcohortlevelp.doc, sideway /// 
replace ctitle(Coef) ///
alpha(0.001, 0.01, 0.05, 0.1) symbol(***, **, *, +)

